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ABSTRACT 

Aims. We present the results of hydrodynamic simulations of the growth and orbital evolution of giant planets embedded in a protoplanetary 
disk with a dead-zone. The aim is to examine to what extent the presence of a dead-zone affects the rates of mass accretion and migration for 
giant planets. 

Methods. We performed 3D numerical simulations using a grid-based hydrodynamics code. In these simulations of laminar, non-magnetised 
disks, the dead-zone is treated as a region where the vertical profile of the viscosity depends on the distance from the equatorial plane. We 
consider dead-zones with vertical sizes, Hdz, ranging from to H^z = 2.3H, where H is the disk scale-height. For all models, the vertically 
integrated viscous stress, and the related mass flux through the disk, have the same value (equivalent to 10"** Mgyr"'), such that the simulations 
test the dependence of planetary mass accretion and migration on the vertical distribution of the viscous stress (and mass flux). For each model, 
an embedded 30 Mg, planet on a fixed circular orbit is allowed to accrete gas from the disk. Once the planet mass becomes equal to that of 
Saturn or Jupiter, we allow the planet orbit to evolve due to gravitational interaction with the disk. 

Results. We find that the time scale over which a protoplanet grows to become a giant planet is essentially independent of the dead-zone 
size, and depends only on the total rate at which the disk viscously supplies material to the planet. For Saturn-mass planets, the migration rate 
depends only weakly on the size of the dead-zone for Hoz < 1.5//, but becomes noticeably slower when Hdz = 2.3//. This effect is apparently 
due to the desaturation of corotation torques which originate from residual material in the partial-gap region. For Jupiter-mass planets, there is a 
clear tendency for the migration to proceed more slowly as the size of the dead-zone increases, with migration rates differing by approximately 
40 % for models with Hoz = and Hoz = 2.3//. 

Conclusions. Our results indicate that for disks models in which the mass accretion rate has a well defined value, the accretion and migration 
rates for Saturn- and Jovian-mass planets are relatively insensitive to the presence and size of a dead-zone. 
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1. Introduction 

' Multiple observations of star-forming regions have revealed 
that young stars are surrounded by circumstellar disks com- 
posed of gas and dust (e.g. Beckwith 1996, Sicilia-Aguilar et 
al. 2006), which are thought to provide the necessary mate- 
rial for planet formation. These protoplanetary disks generally 
show evidence for gas accretion onto the central star, with a 
typical mass flow being ~ 10"** M© yr"' (Sicilia-Aguilar et al. 
2004; Gullbring et al. 1998; Haisch et al. 2001). It has long 
been recognized that such a value for the accretion rate re- 
quires a source of anomalous viscosity, probably originating 
from turbulence, to efficiently transport the angular momentum 
outward (Shakura & Sunyaev 1973). 

Although instabilities such as the baroclinic instability 
(Klahr & Bodenheimer 2003) or the Rossby-wave instability 
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(Lovelace et al. 1999) may result in turbulence, a significant 
number of studies (Hawley et al. 1996; Brandenburg et al. 
1996) have shown that the non-linear outcome of the magne- 
torotational instability (Balbus & Hawley 1991; hereafter MRI) 
is MHD turbulence with an effective viscous stress parameter 
ff ~ 5 X 1 0""*, providing the necessary angular momentum trans- 
port to match the observed accretion rates onto T-Tauri stars 
(Hartmann et al. 1998). 

Although the MRI has been proven robust in generating MHD 
turbulence, the fact that the ionisation fraction is expected to 
be low in protoplanetary disks (Blaes & Balbus 1994) raises a 
number of question about the applicability of the MRI in these 
environments. That is the reason why it has been suggested 
(Gammie 1996) that protoplanetary disks may have both mag- 
netically active zones near the surface, where cosmic ray ion- 
isation enables the MRI to develop, and a "dead-zone" close 
to the disk midplane where the flow remains laminar. This 
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layered-disk picture was confirmed by Fleming & Stone (2003) 
who performed MHD simulations of stratified disks in which 
the magnetic resistivity decreases as a function of height. Over 
the past few years, many studies have focused on the struc- 
ture of the dead-zone by using different chemical reaction net- 
works. Fromang et al. (2002) showed that the size of the dead- 
zone can decrease in the presence of heavy metals (such as 
magnesium) due to charge-transfer reactions. Ilgner & Nelson 
(2006b) examined the role of turbulent mixing, and showed 
that in the presence of magnesium, such a process can enliven 
the dead-zone beyond a few AU. This was confirmed by Ilgner 
& Nelson (2008) who performed multifluid MHD simulations 
and who showed that despite the addition of gas-phase mag- 
nesium, dead-zones in protoplanetary disks typically persist at 
distances R < 5 AU, with potentially important consequences 
for planet formation in these regions. Calculations which ex- 
amine the sizes of dead-zones including gas-grain chemistry, 
show that the presence of a modest population of sub-micron 
sized grains cause the planet forming regions of protoplanetary 
disks to host significant dead-zones (Sano et al. 2000; Ilgner & 
Nelson 2006a; Turner & Drake 2009). 

So far, the effect of a dead-zone on planet formation has re- 
ceived little attention, in large part because of the computa- 
tional expense of performing large-scale 3D MHD simulations 
of turbulent protoplanetary disks. A plausible scenario for the 
formation of planets involves the following steps: i) coagu- 
lation and settling of dust in the disk midplane, followed by 
the growth of km- sized planetesimals; ii) runaway growth of 
planetesimals (e.g. Greenberg et al. 1978; Wetherill & Stewart 
1989) into ~ 10"^ embryos; iii) oligarchic growth of 
these embryos (e.g. Kokubo & Ida 1998, 2000; Leinhardt & 
Richardson 2005) into planetary cores. Planetary cores form- 
ing oligarchically beyond the snow-hne are expected to have 
masses ~ 10 (Thommes et al. 2003) and consequently 
are able to accrete a gaseous envelope to become giant plan- 
ets (Pollack et al. 1996). 

The impact of a dead-zone on dust settling was examined 
by Fromang & Papaloizou (2006) who found a tendency for 
thinner dust sub-disks to form in the presence of a dead-zone. 
Lyra & al. (2008, 2009) found that a Rossby-wave instability 
can be triggered at the border of a dead-zone where the surface 
density is significantly enhanced, creating vortices which are 
efficient in trapping soUds and forming planetary embryos. 
The consequence of a dead-zone on planet migration was 
studied by Matsumura et al. (2003). These authors found 
that a dead-zone significantly slows down the migration of 
low-mass planets undergoing Type I migration. This arises 
because the dead-zone creates a jump in surface density, which 
changed the balance between inner and outer Lindblad torques 
in their model. A surface density transition, however, can also 
create a planet trap where the corotation torque exerted on the 
protoplanet equals the differential Lindblad torque (Masset 
et al. 2006). Concerning gap-opening giant planets which 
migrate on a viscous timescale, Matsumura et al. (2003) found, 
unsurprisingly, that a dead-zone slows down Type II migration 
due to the low value for the viscosity there. However, it should 
be noted that this work did not consider how the existence of a 
Uve-zone near the disk surface affects this latter result. 



In this paper, we present the results of 3D hydrodynamical sim- 
ulations of giant-planets embedded in layered protoplanetary 
disks and in which the dead-zone is modelled using a vertical 
profile for the laminar viscosity which increases as a function 
of height in the disk. The aim of this work is to investigate 
how the evolution depends on the size of the dead-zone Hoz, 
subject to the condition that all our disk models have the same 
vertically integrated viscous stress and associated radial mass 
flux. To address this issue, we consider a model in which a 
30 M® protoplanet is embedded in the disk on a fixed circular 
orbit, and can slowly accrete gas until its mass nip becomes 
the same as that of Saturn (mp = 1 Ms) or Jupiter (mp = 1 
M/). The planet is then released and allowed to evolve under 
the action of disk forces. The results of these simulations 
suggest that for a given accretion rate through the disk, Jupiter 
migrates more slowly as the size of the dead-zone increases. 
However, we find that provided the size of the dead-zone is 
small enough, both the accretion rate onto the planet and the 
migration of Saturn-mass planets depend only weakly on Hoz- 

This paper is organized as follows. In Sect. 2, we describe 
the hydrodynamical model and the numerical setup. In Sect. 3, 
we present the results of our simulations. We finally summarize 
and draw our conclusions in Sect. 4. 



2. Physical model and numerical setup 

In order to study the evolution of a giant planet in a dead-zone, 
which is represented in the simulations as a disk region where 
the viscosity profile depends on the distance from the equato- 
rial plane, we adopt a 3-dimensional disk model. In spherical 
coordinates (r, 6, (p) and in a frame centred on the central star, 
the continuity equation reads: 



dp 

dt 



+ V • (pv) = 



(1) 



where p is the disk density. The equations for the radial, 
meridional, and angular components of the disk velocity v = 
(Vr, ve, v^) are given, respectively, by: 



5(pv,) p(y0 + v^) dp 50) 

^+V.(pV.V)=^ Y^-P^^fr, 



(2) 



djpve) 
dt 

and 



+V-(pv0\) = -^+fe (3) 

r r r or r or 



dt 



+ V • (pv^v) = ■ 



pVrV^ pvev^cote 



1 dp 
r sin 6 dr 



rsinff dr 



(4) 
(5) 



In the above equations, p is the pressure, fa and are re- 
spectively the radial, meridional and azimuthal components of 
the viscous force per unit volume. Expressions for fr, fg and 
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can be found for example in Klahr et al.(1999). O is the gravi- 
tational potential and can be written as: 



CM, 



Gm„ 



+ 0, 



(6) 



where M* and irip are the masses of the star and planet, respec- 
tively, and where e is a softening length. O,,,^ is an indirect term 
arising from the fact that the star-centered frame is not inertial. 
This term reads: 



Gm„ 



ind 



-r ■ Tp -H G 



r dmir') , 



(7) 



where the integral is performed over the volume of the disk. 
2.1. Planet orbital evolution 

In this work the planet can experience the gravitational acceler- 
ation arising from both the disk and the central star. Therefore, 
the equation of motion for the planet is given by: 



dfi 



+ fdp - V(5„ 



where fdp is the force due to the disk which is defined by: 

p(r)(rp - r) 



Idp 



Jv (Ir-TpF 



rdV 



(8) 



(9) 



|2 + ^2)3/2' 

Note that we exclude the material contained in the planet 
Hill sphere when calculating the gravitational force acting on 
the planet. Moreover, in the simulations presented here, the 
smoothing length e is set to the diagonal length of one grid 
cell. 

3. Numerical setup 

3.1. Numerical method 

The 3D hydrodynamical simulations presented in this paper 
have been performed using the NIRVANA code (Ziegler & 
Yorke 1997), which is basically a grid-based code which com- 
putes spatial derivatives using finite differences. The advec- 
tion scheme is based on the monotonic transport algorithm and 
since it employs a staggered mesh, the numerical method used 
in this code is spatially second-order accurate. Further details 
about NIRVANA can be found for example in De Val-Borro et 
al. (2006). 

We adopt computational units in which the mass of the central 
star is - 1, the gravitational constant is G = 1 and the 
cylindrical radius R - \ corresponds to the initial orbital radius 
of the planet ao- In the following, we report our results in units 
of the initial orbital period of the planet P - 27rf2o', where 

Qo = ^GMJal 

For most of the calculations, we employ A^^ = 128 radial grid 
cells uniformly distributed between r,„ = 0.4 to rout - 2.5, 
- 384 grid cells in azimuth and Ne - 64 meridional grid 
cells with 6 lying in the range 80° to 90°. Some low resolution 
runs using Nr = 64, A^^ - 192 and Ng = 64 grid cells have 
also been performed. These are long-term simulations, aimed 
at giving the embedded giant planet sufficient time to create a 
gap in the disk. 




Fig. 1. This figure shows the a-parameter as a function of the 
distance from the equatorial plane for the different sizes of the 
dead-zone we considered, namely for Hoz/H = (black line), 
Hdz/H = 0.5 (redhne), Hoz/H = 1 (green Une), Hdz/H = 1.5 
(blue line) and Hoz/H = 2.3 (orange line). 

3.2. Initial conditions 

In the disk model that we adopt for the simulations presented 
here, the aspect ratio is constant and set to h - H/R = 0.05, 
which is a typical value for protoplanetary disks. The initial 
density profile is chosen such that the vertical stratification of 
the disk follows the condition for hydrostatic equilibrium and 
is defined by: 

where R is the cylindrical disk radius, Z is the distance from 
the midplane and where po and Rq are constants. Consequently, 
the disk surface density 1. - J pdz can be written as S = 

I.o(Ro/Ry^'^, with So = y/lnhRopQ. Here, we set Rq = 1 and 
Po is chosen such that the disk contains ~ 3 x lO^^M© within 
the computational domain, which corresponds to ~ 0.02 Mq 
interior to 40 AU. 

The anomalous viscous stress arising from MHD turbulence is 
parameterized using the standard 'alpha' prescription for the 
disk viscosity v - ac^H (Shakura & Sunyaev 1973), where c, 
is the local speed of sound . From the results of previous MHD 
simulations (Papaloizou & Nelson 2003; Fromang & Nelson 
2006), such a prescription is expected to provide a reasonable 
description of the mean flow in the disk. It is now understood 
that low mass planets are subject to stochastic forcing when 
embedded in turbulent disks (e.g. Nelson & Papaloizou 2004; 
Nelson 2005), and so the adoption of an anomalous viscous 
stress to model the effects of turbulence may seem inappropri- 
ate at first sight. But, in this work, we are mainly interested in 
the migration and growth of more massive, gap forming plan- 
ets where the primary role of turbulence is to supply gas to 
the vicinity of the planet through an accretion flow generated 
by the turbulent stresses, and the local interaction between the 
planet and turbulent density fluctuations is less important (e.g. 
Nelson & Papaloizou 2003; Papaloizou, Nelson & Snellgrove 
2004). 
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Although it remains to be proven that for a disk with a dead- 
zone, the stresses arising from MHD turbulence behave like 
a locally varying laminar viscosity, a is chosen here to be Z- 
dependent and such that it can vary from very small values in 
the disk midplane to much higher values near the disk surface. 
More precisely, in the calculations we performed, the vertical 
profile adopted for a was as follows; 

a(Z) = - Uoaz - aoz) tanh ^ A~^) ^ ^ "'^^] ' '"^^^ 

where Hoz is the size of the dead-zone and where oaz (resp. 
aoz) is the value for the a-parameter in the active (resp. 
dead) region. In the previous equation, A denotes the size 
of the transition from the dead to the active zone and is set 
to A = 0.08// in the work presented here. In order to study 
the influence of the dead-zone size on the evolution of giant 
planets, we have considered different values for Hoz lying 
in the range to 2.3 H. In the case where Hoz = 0, a is 
constant in the entire disk and is set to ao = 1.4 X 10"^, which 
corresponds to a constant accretion rate through the disk 
Mo = 3naoI,oh^. For other values of Hoz, we adopt a similar 
value for the vertically averaged mass flow through the layered 
disk M - 3n jv(z) pdz = Mo and we calculate the value for 
Oaz accordingly, assuming that the residual viscosity in the 
dead-zone is aoz = 1.0 x 10"^. For the different values of Hoz 
we considered, the coiTesponding value for uaz can be found 
in Table 1 and Fig. [T] displays the corresponding a-parameter 
as a function of the distance from the disk equatorial plane . 

The planet is initially placed on a circular orbit at a = 1 , and 
our main interest is to examine the orbital evolution of planets 
with masses nip - 3 x 10""* and lO""' M,, corresponding to 
Saturn-mass (Ms) and Jupiter-mass (My) planets, respectively. 
In order to conserve mass as the planet grows and opens a gap, 
we proceed in two steps: 

i) we consider a 30 M^ protoplanet held at a = 1 on a fixed 
circular orbit. This body is allowed to accrete gas from the disk 
on a timescale coiTesponding to ~ 1500 orbits. Due to the very 
long run times required for the planet mass to reach 1 My, we 
have performed, for this step, low-resolution simulations using 
Nr = 64, A^^ = 192 and Ne = 64 grid cells. 

ii) once the planet mass has reached 1 Mg or 1 My, we 
restart the calculations, but with the planet being allowed to 
migrate due to its interaction with the disk, and accretion being 
switched off. For this second step, we have performed higher 
resolution simulations using A^,. = 128, - 384 and Ne - 64 
grid cells. The new values for the disk physical quantities were 
computed from the old ones using bi-linear interpolation. 

3.3. Boundary conditions 

At the inner edge of the computational domain, we model ac- 
cretion onto the central star by using a 'viscous' boundary con- 
dition, for which we set the radial velocity in the innermost 
cells to Vr = /3vyi,c(n„), where v,.,jc('-m) = -3v/2r,„ is the typ- 
ical disk inward drift velocity due to viscous diffusion, and j6 
is a free parameter which is set to y6 = 5 in this work (Pierens 



Table 1. Dead-zone parameters for the different models con- 
sidered. 







Cl-DZ 
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Fig. 2. This figure shows the evolution of the planet mass as 
a function of time for the different sizes of the dead-zone we 
consider, namely for HozlH - (black line), HozlH - 0.5 
(red line), HozlH - 1 (green line), HozlH - 1.5 (blue line) 
and HozlH - 2.3 (orange line). 

& Nelson 2008). At the outer edge of the computational do- 
main, and at the lower meridional boundary as well, we prevent 
mass loss from the disk by employing reflecting boundary con- 
ditions. Assuming that the disk is symmetric with respect to its 
midplane, we also impose a symmetry boundary condition at 
the upper meridional boundary, which corresponds to the disk 
midplane. Consequently, we consider only the upper half of the 
disk in the simulations described below. 

4. Results 

4.1. Effect of the presence of a dead-zone on 
planetary growth 

For each value of Hoz we consider, we have performed a 
simulation in which a 30 Mg planet held on a circular orbit can 
accrete gas from the disk until its mass reaches 1 My . In the 
simulations presented here accretion is modelled by removing 
at each time-step a fraction of the gas located inside the Roche 
lobe of the planet and then adding the corresponding amount 
of matter to the mass of the planet (e.g. Kley 1999; Nelson 
et al. 2000). The e-folding time for gas accretion was chosen 
to be the orbital period of the planet, which corresponds to 
the maximum rate at which the planet can accrete gas (Kley 
1999). The choice of a 30 M® protoplanet to act as a seed onto 
which gas can accrete is broadly consistent with evolutionary 
models of gas giant planets forming in protoplanetary disks, 
which suggest that rapid gas accretion occurs once the planet 
mass reaches 30 - 40 M® (Papaloizou & Nelson 2005). 
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The time evolution of the planet mass for each model is dis- 
played in Fig. 12] We see that the evolution of the planet mass 
depends only weakly on the size of the dead-zone and proceeds 
similarly in each case. The early evolution involves rapid gas 
accretion with the planet growing to become a Saturn-mass 
planet in ~ 100 orbits. As the planet mass increases, non-linear 
effects can tidally truncate the disk, leading to a decrease in the 
accretion rate such that it takes ~ 1500 orbits for a Saturn-mass 
planet to reach a Jovian mass. 

It is worth noting, however, that although at earlier times the 
accretion rate should not depend strongly on the size of the 
dead-zone, since the planet is able to accrete material which is 
on neighbouring orbits, we expect this to be not true from the 
time the gravitational torque is largely deposited in the disk in 
the near-vicinity of the planet. This occurs when the thermal 
criterion for gap opening is satisfied (Lin & Papaloizou 1986; 
Crida et al. 2006), namely when the radius of the Hill sphere 
exceeds the disk semi-thickness: 



H Dz /H=2.3, m 



3M, 



1/3 



> H. 



(12) 



In the simulations presented here, this happens roughly when 
nip > 3 X 10 From this time, the planet may be able to open 
a gap in the disk provided that the gravitational torques over- 
whelm the viscous forces (Bryden et al. 1999; Papaloizou et al. 
2006), a condition which can be expressed as: 



nip 40v 
— ^ > — , 



(13) 



where Q^, is the angular velocity of the planet. Since the value 
for the viscosity in the active region is correlated with the size 
of the dead-zone, the previous equation shows that the gap 
structure and consequently the accretion rate onto the planet 
may depend strongly on Hdz for irip > 3 x 10 Indeed, 
considering for instance the case nip = 1 Mj, Eg. [T3] predicts 
that gap opening should occur for a < 0.01. Therefore, in the 
cases where Hoz/H = 1.5 and HqzIH - 2.3, we expect a 
I Mj planet to be able to open a gap only in the dead-zone, 
possibly leaving the gas flowing over the top of the planet in 
the active region if the accretion flow there is rapid enough. 
Snapshots of the disk density in the equatorial plane, and at the 
disk surface, are depicted in Fig. [3] These plots correspond to 
a time when the mass of the accreting protoplanet has reached 
1 Mj and are for a model in which Haz - 2.3//. These reveal 
that the planet truncates the disk throughout its vertical extent, 
despite the high value for the viscosity in the active region. 
As can be seen in Fig. 21 such a process arises because the 
planet is able to pull some of the gas from the live-zone down 
in toward the midplane as the latter tries to flow over the top 
of the planet. This suggests that both the gap structure and the 
accretion rate onto the planet do not depend strongly on the 
presence of a dead-zone, but rather on the vertically-integrated 
accretion rate of gas through the disk. In that case, we would 
expect the results from these 3D runs to not significantly differ 
from 2D calculations for a given value of the accretion rate M. 




Fig. 4. This figure shows, for HozlH - 2.3, details of the flow 
in the neighbourhood of the accreting planet when its mass has 
reached 1 Mj. The velocity field is represented in the {R,z) 
plane and for 0.6 < /? < 1.5. 



To investigate this issue in more detail, we have performed 
a suite of 2D runs with varying values of a from ff = to 
a - 1.4 X 10"^ and for which the initial surface density at 
the planet position is 2(/?o) = Sq. An additional 2D calcula- 
tion with a - 2.8 x lO""* and S(/?o) = 0.52o has also been 
performed. The results of these calculations are shown in Fig. 
|5] Also displayed are the results of two 3D simulations corre- 
sponding to: i) a fully active disk with a - 1.4 x 10"^ through- 
out; ii) a dead disk with a - Q everywhere. As expected, 2D 
runs with a i^O predict that the protoplanet grows faster as the 
viscosity increases. The run with a - 2.8 x 10"^ shows slower 
growth at earlier times since the surface density is smaller in 
that case. However, once the planet opens a gap in the disk, it 
is clear that the growth rate becomes similar to that obtained 
from the run with a - 1.4 x 10"^ and for which the value for 
the mass accretion rate through the disk is the same. In the case 
where a = 0, viscous supply of gas through the disk is inhib- 
ited and the planet mass can saturate once the feeding zone is 
empty of gas, which occurs when nip ~ Q.SMj for the disk 
parameters used in this work. Clearly, for both a = and 
a - 1.4 X 10"\ there is a good agreement between 2D and 
3D calculations. This confirms that the accretion rate onto a 
gap-opening planet does not depend strongly on the details of 
the flow in the (R, z) plane, but only on the vertically-averaged 
accretion rate through the disk. 



4.2. Effect of a dead-zone on the orbital evolution of 
Saturn-mass planets 

In order to study the evolution of a Saturn-mass planet embed- 
ded in a layered protoplanetary disk, we restarted the simula- 
tions described in the previous section when the mass of the 
accreting planet had reached the relevant value. We then re- 
leased the planet and let it evolve freely under the action of 
disk torques. 

The evolution of Saturn's semi-major axis for the different 
models is shown in the upper panel of Fig.|6l For each value of 
HozlH, Saturn does not undergo runaway migration (Masset 
& Papaloizou 2003) because the disk is not massive enough. 
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Fig. 3. Left panel: this figure shows a snapshot of the disk density in the equatorial plane for Hdz/H - 2.3. Here, a planet held 
on a circular orbit at a = 1 and with initial mass nip - 30 M© has grown to become a Jupiter-mass planet. Right panel: same but 
at the disk upper surface. 
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Fig. 5. This figure shows the evolution of the planet mass as a 
function of time for 2D simulations in which the (T-parameter 
was varied from a = to a - 2.8 x lO"-'. Diamonds display 
for comparison the results of a 3D run with Hdz - and a = 
1 .4x10"^' whereas triangles display the results of a 3D run with 
a = 0. 

Instead, the migration rate is intermediate between the Type II 
and Type I regimes. This arises because the mass of Saturn 
is insufficient to enable a clean gap to form in the disk, as 
illustrated in Fig. [8] This figure displays, for the model with 
Hdz/H = 2.3, a snapshot of the Saturn-induced gap structure 
in both the equatorial plane and at the disk surface. It is worth- 
while noticing that gap opening occurs in the live-zone, de- 
spite the high value for the viscosity there, indicating that both 
Jupiter and Saturn can pull gas from the active region down to- 
ward the disk midplane. 

For models with Hoz /H < 1.5, the evolution of the system pro- 
ceeds quite similarly, with a slight tendency for the migration 
rate to decrease as the size of the dead-zone increases. For the 
calculation with Hdz/H = 2.3 however, we see that in ~ 200 
orbits the semi-major axis has decreased by an amount that is 
nearly 20% less than for the other cases. This can be confirmed 



by examining the total disk torques exerted on the planet, which 
are displayed for each model in the lower panel of Fig. |6] It is 
well known that the total tidal torque, F, exerted on the planet 
can be written as the sum of a differential Lindblad torque, F^, 
and a total corotation torque, Fc (Tanaka & Ward 2002). In con- 
trast to Fi, which is basically independent of viscosity (Meyer- 
Vernet & Sicardy 1987; Papaloizou & Lin 1984), Fc depends 
strongly on v and can saturate in the low-viscosity limit. In the 
present case, Fc can be decomposed as: 



F?2+F^^ 



(14) 



F^. ) is the corotation torque due to the li- 



where F^^ (resp. 
brating fluid elements which originate from the dead (resp. ac- 
tive) region. In the dead-zone, F^^ is expected to saturate on 
a timescale corresponding to the outermost horseshoe turnover 
time (Masset 2002) 



3jc, p 



(15) 



In the previous equation, X j is the horseshoe zone width which 
can be approximated by Xs = 2A5a(q/3y^^ (Masset et al. 
2006), where q is the planet to central star mass ratio. Here, 
this gives ths ~ 71 Q^', which is consistent with the timescale 
over which the total torques reach a constant value in Fig. |6] 
In the active region, however, we would expect the corota- 
tion torque F^^ there to saturate only if the viscous time-scale 
T,, = jc, /3v across the horseshoe region exceeds ths ■ Assuming 
that Xs is constant over the vertical extent of the disk, and for 
the disk parameters used in this work, this suggests that the 
corotation torque saturates in the active region provided that 
fffl < 2.5 X 10"^. Therefore, we would expect the corotation 
torque to saturate in models with Hdz/H < 1.5. For the cal- 
culation with Hdz/H = 2.3 however, we would expect the 
corotation torque originating from the live-zone to be positive 
since for our disk model, the disk surface density decreases 
/? Such a scenario is consistent with the fact that the to- 
tal torques exerted on Saturn do not depend on the size of the 
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Fig. 6. Upper panel: this figure shows the evolution of Saturn's 
semi-major axis for the different sizes of the dead-zone we con- 
sider, namely for HqzIH - (black line), HazlH - 0.5 (red 
line), HazlH - 1 (green line), Hoz/H =1.5 (blue line) and 
Hdz/H = 2.3 (orange line). Lower panel: this figure shows the 
corresponding total torque exerted on the planet. 



dead-zone for Hoz/H - 1.5, and are weaker in the case where 
Hdz/H = 2.3 (see Fig. |6]l. For the run with Hpz = 2.3//, 
the computed mass of gas material located inside the coro- 
tation region is ~ 8 x 10"^ in our units with about 3 % of 
this mass located in the active zone, which appears to be suffi- 
cient to provide a non negligible contribution to the total coro- 
tation torque exerted on the planet. In order to clearly demon- 
strate that the difference in the torques arise from the viscously- 
evolving corotation region, we show in the upper panel of Fig. 
|2]the radial torque distribution r(r) for the runs with Hdz - H 
and Hoz = 2.3// and for a planet located at a = 0.96. In com- 
parison with the simulation in which Hdz = H, the run with 
Hdz = 2.3// exhibits a higher positive bump at r ~ 0.89 
and a smaller negative bump at r ~ 1.04, which therefore 
corresponds to a slightly smaller differential Lindblad torque. 
For the run with Hdz = 2.3//, the positive torque excess for 
0.9 < r < 0.94 can be explained by looking at the vertical 
torque distribution r(z) exerted by the disk region located be- 
tweena-Xs anda+x, and whichis displayed in the lower panel 
of Fig. |7] For this simulation, the torque becomes positive be- 
yond Z > 1 .7//, which is clearly not the case for the run with 
Hdz - H. This unambiguously confirms that the corotation 
torque exerted by the viscously-evolving layers is unsaturated 
in the run with Hdz = 2.3// and tends to slow down migration. 
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Fig. 7. Upper panel: this figure shows the radial distribution 
of the torque for the runs with Hdz - H and Hdz/H = 2.3. 
The planet is here located at a = 0.96. Lower panel: this figure 
shows the vertical distribution of the torque in the upper disk 
layers. 



4.3. Orbital evolution of Jupiter-mass planets 

We now turn to the question of how the orbit of a Jupiter-mass 
planet evolves in a dead-zone. Here again, to address this is- 
sue, we restarted the simulations presented in Sect. 14.11 once 
the mass of the accreting protoplanet has reached nip - \ Mj, 
but we now let the planet evolve under the influence of the disk 
for ~ 500 orbits. 

The time evolution of Jupiter's semi-major axis is shown in the 
upper panel of Fig. |9] As expected, the planet migrates on a 
timescale coiTesponding to Type II migration in the case where 
Hdz/H - 0. For models with Hdz + however, we see here 
that, compared with the Saturn case, the dependency of the mi- 
gration rate upon the value for Hdz is stronger, with a clear 
tendency for the semi-major axis to decrease more slowly as 
the size of the dead-zone increases. Examination of the total 
disk torques, which are represented in the lower panel of Fig. 
|9l clearly reveals that these decrease in magnitude with increas- 
ing the size of the dead-zone. The outer and inner torques ex- 
erted on the planet are displayed in the upper and lower panels 
of Fig. [TOl respectively. Interestingly, and despite a weak de- 
pendency on the value for Hdz, the presence of a dead-zone 
reduces the effects of the outer torques exerted on the planet. 
As illustrated in Fig.[TT] which shows the disk surface density 
as a function of radius (at an azimuthal position corresponding 
to that of the planet), this arises because the surface density at 
the outer edge of the gap is greater in the case where Hdz - 



A. Pierens and R.P. Nelson: On the growth and orbital evolution of giant planets in layered protoplanetary disks 



H /H=2.3 , Disk equatorial plane, mp=3. 1 ' 



Hd,/H=2.3 , Disk surface, m„=3.10 '* 




l.1e+00 > 




500 



Fig. 8. Left panel: this figure shows,for Hdz/H = 2.3, a snapshot of the disk density in the equatorial plane at a time where a 
planet held on a circular orbit at R - I and with initial mass nip - 30 has grown to become a Saturn-mass planet. Right 
panel: same but at the disk surface. 



torque exerted on the planet in this case (see the lower panel 
of Fig. [Tol l. It should be noted that for models with Hdz + 0, 
and provided that the viscosity is high enough in the live-zone, 
a migration rate which is lower than the inward drift velocity 
in the active region can make the disk material from this re- 
gion flow across the gap. This not only enables the inner disk 
to be continuously supplied with gas but also creates a positive 
corotation torque on the planet as the gas flows through the gap 
(Masset & SneUgrove 2001; Morbidelh & Crida 2007). 

4.4. Planet evolution in an alternative dead-zone 
model 

The prescription that we have used to model the variation of 
viscous stress in disks with dead-zones depends only on the 
height of the fluid from the disk midplane. It is possible that 
some of our results may be affected by this prescription, since 
the spatial extent of a dead-zone should probably depend on the 
local column density (assuming that external sources of ionisa- 
tion such as cosmic rays or X-rays are responsible for ionising 
the disk). The assumption of a stationary profile for a is rea- 
sonable in a disk which remains largely unperturbed, but may 
become inaccurate when gaps form due to the presence of plan- 
ets, since material near the midplane located in the gap should 
become magnetically active. To further investigate this issue, 
we have performed additional simulations using the following 
vertical profile for a: 
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Fig. 9. Upper panel: this figure shows the evolution of Jupiter's 
semi-major axis for the different sizes of the dead-zone we con- 
sider, namely for Hoz/H - (black line), Hoz/H = 0.5 (red 
line), Hdz/H - 1 (green line), Hoz/H =1.5 (blue line) and 
Hdz/H - 2.3 (orange line). Lower panel: this figure shows the 
corresponding total torque exerted on Jupiter. 

compared to cases where Hdz > 0, increasing the magnitude 
of the outer torques. Moreover, Fig.[TT]reveals that the density 
in the inner disk is much more higher for the simulation with 
Hdz - 2.3//, which consequently favours the higher positive 



a - ao(R) exp 



(16) 



where r/iz) - pdz' is the column density at height z, and 
where ao is determined in such a way that the vertically aver- 
aged mass flow through the disk is Mq at each radius. In the 
previous equation, we set rjc = 8.35 x 10"^ which coiTesponds 
to the initial value for the column density at altitude z = 1 .5H. 
The upper panel of Fig. [12] shows the azimutally-averaged 
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Fig. 10. Upper panel: this figure shows the time evolution of 
the outer torque exerted on Jupiter, for the different sizes of 
the dead-zone we consider,namely for HqzIH - (black line), 
HdzIH = 0.5 (red line), HdzIH = 1 (green line), HdzIH = 1.5 
(blue line) and HdzIH - 2.3 (orange line). Lower panel: same 
but for the inner torque. 



vertical profile of the a-parameter aiR - \ for a simulation in 
which an initial 30 protoplanet grows to become a Jupiter 
mass planet, and for which a is given by the previous equation. 
The solid line depicts the initial profile for a, whereas the 
dashed and dot-dashed lines correspond to times at which 
nip - 3 y. lO""* and nip - 10"^, respectively. As the planet 
grows and is able to open a gap in the disk, the vertical profile 
of a at the position of the planet relaxes toward a uniform 
distribution with a - 1.4 x lO^^* throughout the vertical extent 
of the disk. Although not presented here, it is worthwhile 
noting that for such a simulation, the time evolution of the 
planet mass is almost indistinguishable from that obtained 
using the prescription of Eq. [TT] as the value for the mass 
accretion rate is the same. 

For this model, the middle and lower panels of Fig. [12] 
respectively, display the evolution of the semi-major axes 
for the Saturn and Jupiter-mass planets. Compared with the 
previous model with Hqz - 15//, the migration of Saturn 
proceeds very similarly. This arises because in both cases, 
the value for the laminar viscosity in the active region is not 
high enough for the corotation torque to be unsaturated. In the 
Jupiter case, the difference between the two models after 200 
orbits is again marginal despite the fact that for the new model 
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Fig. 11. This figure shows, for the simulations in which nip — 
10"^, the disk surface density as a function of radius and for an 
azimuthal position corresponding to that of the planet. Different 
sizes of the dead-zone are considered, namely HdzIH - 
(black line), HdzIH = 0.5 (red line), HdzIH = 1 (green line), 
HdzIH -1.5 (blue line) and HdzIH - 2.3 (orange line). 



a is vertically constant at the position of Jupiter (see the upper 
panel of Fig. [T2l i. This supports the idea that the prescription 
we adopt for a has little influence on the subsequent evolution 
of the planet. 



5. Discussion and conclusion 

In this paper we have presented the results of hydrodynamic 
simulations aimed at studying the effect of a dead-zone on the 
formation and evolution of giant planets. For layered disks we 
model the dead-zone as a region where the viscosity profile 
depends on the distance from the equatorial plane, and we 
considered different values for the dead-zone size, ranging 
from Hdz = to Hdz - 2.3H. For each value of Hdz, we 
focused on a model in which a 30 protoplanet accretes 
gas from the disk until its mass reaches either nip = 1 Ms or 
nip - 1 Mj. Once the final planet mass has been attained, we 
assumed that accretion stops, and let the planet evolve under 
the action of disk torques. 

The results of our simulations indicate that, for a given mass 
accretion rate through the disk, the timescale over which a 
30 Me protoplanet grows to become a giant planet does not 
depend on the size of the dead-zone. This occurs because 
the planet is able to pull some of the gas from the live-zone 
down toward the midplane, from where it can be accreted, 
as the latter flows into the gap region. However, we find that 
the migration of giant planets can be modestly affected by 
the presence of a dead-zone. Indeed, there is a clear tendency 
for Jupiter-mass planets to migrate more slowly as the size 
of the dead-zone increases. For Saturn-mass planets, we find 
that although migration depends only weakly on the size 
of the dead-zone for models with Hdz ^ 1.5//, it proceeds 
more slowly in the case where Hdz = 2.3//. This is due to 
the fact that in the latter case, the viscosity in the live-zone 
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Fig. 12. Upper panel: this figure displays the (T-parameter 
azimutally-averaged vertical profile alR - \ for a simulation in 
which an initial 30 Me protoplanet grows to become a Jupiter- 
mass planet and in which a is given by Eq. [16] Different times 
are considered corresponding to times at which nip - 9.10"^ 
(solid line), nip - 3 x 10""* (dashed line) and nip = 10"^ (dash- 
dotted line). Middle panel: this figure shows the time evolution 
of Saturn's semi-major axis, for a model in which a is given by 
Eq. nn with Hoz/H -1.5 (blue line) and for a simulation in 
which a is given by Eq. [T6l (dashed line). Lower panel: same 
but for Jupiter. 

is high enough for the corotation torque to become unsaturated. 

Although the present study suggests that a dead-zone does 
not have a significant impact on the growth and migration of 
giant planets, it will be of interest to examine this issue using 
three-dimensional MHD simulations, which include a dead- 
zone, in order to confirm (or refute) our findings. Of particular 
interest will be the issue of what happens when gap opening al- 
lows disk material in the vicinity of the planet to become mag- 
netically active. We will address this issue in a future paper 
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